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Abstract 

The non-LTE state of the hydrogen gas in isothermal transonic decretion disks 
around a BIV star has been calculated by an iterative method in order to explore 
the basic physical process in the disk. This dynamical model is characterized by 
a density law in the equatorial plane of p{R) oc R~^'^. The continuous radiation is 
calculated with the A iteration in the integral form, while we adopt a single-flight 
escape probability for lines. We describe the non-LTE state, the radiation flow and 
conversion in the disk. We conclude that the stellar Balmer continuum plays a key 
role in the non-LTE state of the disk. The examination of the local energy gain and 
loss suggests that the disk temperature has double minima along the equatorial plane 
in the optically thick case: the intermediate region caused by deficient ultraviolet 
radiation and the outer Lyman a cooling region. We have also calculated some 
observable quantities, such as the spectral energy distribution, the UBV colors, the 
infrared excess and the Balmer line profiles. Our calculations with the mass loss rate 
less than 10~^^Mq yr~^ reproduce the observed continuum quantities. However, we 
could not get large Ha emission strength observed in Be stars. We suggest that the 
density gradient of the Be star disk is slower than that of the isothermal decretion 
disk. 

Key words: stars: atmospheres — stars: circumstellar matter — stars: emission- 
line. Be — line: profiles — radiative transfer 

1. Introduction 

Be stars are classified as main sequence or giant B stars which have ever shown the 
emission in Balmer lines. Struve (1931), in his pioneering work on Be stars, proposed flattened 
equatorial disks around Be stars from the correlation between the emission line widths and the 
projected rotational velocities of underlying stars. Coyne & Kruszewski (1969) supported this 
view from the polarimetric observations. While, the rotation law in the disk had been unsettled 
for a long time. It seems that it has been recently established that the disk is rotating with a 
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Keplerian law, because the one-armed oscillation, which is excited only in the Keplerian or the 
near Keplerian disk, can explain the observational facts such as the V/R variation (Okazaki 
1997) and the Balmer progression (Hirata & Okazaki 2000). The most confident evidence 
for the Keplerian disk is the interferometric detection of a size difference of the Ha emission 
region in the anti-direction measured from the central star, which correlates with the V/R 
variation (Vakili et al. 1998; Berio et al. 1999). Hummel & Vrancken (2000) obtained the 
same conclusion from the comparison between the observed and calculated emission profiles. 
The decretion disk, which was introduced first by Lee, Saio & Osaki (1991) and later developed 
by Porter (1999) and Okazaki (2001), gives the physical basis of the Keplerian disk, though 
the mechanism which supplies the mass and angular momentum is not yet known. The Be 
star disk model with magnetic field has been devised by Cassinelli et al. (2002), Maheswaran 
(2003) and ud Doula & Owocki (2002). However, Owocki & ud Doula (2003) concluded that 
the simulation with the stellar rotation does not form a magnetically torqued disk with enough 
density. 

It is well known that the Be star disk is highly in non-LTE state, because the density is 
not so high that the coUisional process is effective. It is necessary to solve non-LTE problem 
to derive physical quantities from observational materials. Then, it is highly desirable to 
calculate non-LTE state for the decretion disk model. Several non-LTE calculations have been 
published for the flattened disks with different velocity structures, namely, different density 
structures. Marlborough (1969) and Poeckert & Marlborough (1978) adopted the weighted 
mean of some extreme cases as the non-LTE solution. Millar & Marlborough (1998) and Millar 
& Marlborough (1999) expanded the model to the radiative equilibrium model. The latter took 
into account the diffuse radiation in an approximate way. The model has been developed further 
by Jones, Sigut & Marlborough (2004) and Sigut & Jones (2007), in which they examined the 
metal effect to the thermal structure of the disk. Sigut & Jones (2007) adopted a single-flight 
escape probability in the static plane-parallel geometry for lines and on-the-spot approximation 
for continuous radiation. They also took into account the diffuse continuous radiation in an 
approximate way. Non-LTE calculations were also performed by Kfi'z (1979), Stee et al. (1995), 
and Stee & Bitter (2001), based on Sobolev approximation. This approximation is an adequate 
one when the velocity gradient is large (Hamann 1981; de Koter, Schmutz & Lamers 1993). 
However, the expansion velocity in the decretion disk is very small (Okazaki 1997). Hummel 
(1994) and Hummel & Hanuschik (1997) adopted the equivalent two level approximation and 
calculated Ha line proflles, using A iteration. Carciofl & Bjorkman (2006) performed the 
non-LTE calculation for the decretion disk under the radiative equilibrium by using Monte 
Carlo method. Carciofl et al. (2007) improved their model to take into account the vertical 
temperature structure. They showed that there is a temperature minimum in the several stellar 
radii on the equatorial plane. 

In this paper, we present our results based on the full treatment of the diffuse continuous 
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radiation, as an extension of Iwamatsu & Hirata (2007). In section 2, we describe the central 
star and the disk structure for which we calculate the non-LTE state. In section 3, we describe 
our method of the non-LTE calculation. Section 4 gives the non-LTE state in our model, the 
local energy balance analysis, and also the total energy gain and loss in the whole disk. In 
section 5, we describe the continuum properties such as the spectral energy distributions, the 
UBV colors, and the infrared excess. Section 6 treats the line properties such as the Balmer line 
profiles, the equivalent widths, and the Balmer decrement. The summary and the concluding 
remarks are given in section 7. 



2. Central Star and Decretion Disk 

The spectral type of the central star is assumed to be BIV and its basic parameters 
are chosen as the effective temperature T^s = 24000 K, the radius i?* = 6.9 Rq, and the mass 
= II.OMq, referred to vr Aqr (Bjorkman et al. 2000). The stellar flux is taken from ATLAS9 
(Kurucz 1993) with the other parameters: the gravitational acceleration at the stellar surface, 
logg = 4.0, the solar abundance, and the turbulent velocity of 2 kms~^. For easier integration 
over the frequency, we adopt the stellar flux without line blanketting. Although the effective 
temperature for the stellar flux without line blanketting is actually 26400 K, we expect that the 
influence of line blanketting is small, because the flux difference near each shortward continuum 
limit is small. We assume the central star is spherical without limb darkening. 

As a dynamical model, we adopt an isothermal transonic decretion disk in steady state 
with a viscosity parameter a = 0.1 (Okazaki 2001). This one-dimensional model is character- 
ized by almost Keplerian rotation and small expansion velocity. We assume that the disk is 
axisymmetric about the rotational axis and symmetric about the equatorial plane, and is in 
hydrostatic equilibrium in the vertical direction. It is convenient to adopt a cylindrical coordi- 
nate, whose origin and the polar axis are the center of the star and the rotation axis of the star, 
respectively. Let us denote the distance from the rotation axis by w and the vertical distance 
from the equatorial plane by z in units of the stellar radius. We assume that the expansion 
and rotational velocities at any point are the same as those at w in the equatorial plane. The 
number density N{w,z) is related to the number density along the equatorial plane N{w,0) 
(Marlborough 1969) as 

iV(...)=/V(..0)exp{-i(i-^=^)). (1) 

with 

where /i is the mean molecular weight, which is fixed to 0.5 in our calculations. M^, is the 
stellar mass and Tdisk is the temperature of the isothermal disk. The other notations have 
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Table 1. Physical quantities of our dynamical model sequence. We show the mass loss rate M, the total mass 
of the disk M^^^]^, the disk temperature T^jigij, and the density 0), at ro = 1 and 1000. a{b) means a X lO''. 



model A/ A/disk Tdisk p{^,0) [g/cm ^] 

[Moyr^i] [Mq] [K] zu = 1 m ^ 10^ 

1 1.0(-12) 1.02(-10) 1.6(4) 1.75(-12) 8.60(-24) 

2 2.0(-12) 2.04(-10) 1.6(4) 3.50(-12) 1.72(-23) 

3 5.0(-12) 5.11(-10) 1.6(4) 8.74(-12) 4.30(-23) 

4 l.O(-ll) 1.02(-09) 1.6(4) 1.75(-11) 8.60(-23) 

5 2.0(-ll) 2.04(-09) 1.6(4) 3.50(-ll) 1.72(-22) 

6 5.0(-ll) 5.11(-09) 1.6(4) 8.74(-ll) 4.30(-22) 

7 I.O(-IO) 1.02(-08) 1.6(4) 1.75(-10) 8.60(-22) 

8 2.0(-10) 2.04(-08) 1.6(4) 3.50(-10) 1.72(-21) 

9 I.O(-IO) 1.68(-08) 1.2(4) 3.15(-10) 1.35(-21) 



usual meanings. The number density N{u7,0) is obtained from the column density of Okazaki 
(2001 )'s model under the assumption that the disk is geometrically thin. This assumption is 
given by 

2 
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This dynamical model is characterized by N{to,0) oc The number density A^(ci7,2;) is 

derived from equation (1). Table 1 shows the mass loss rate, M, the total mass of the disk 
-^disk, the assumed disk temperature Tjisk, and the density p(ro,0) at w = l and 1000 in our nine 
models ^ The number density A^(1,0) is given by A^(1,0) = p(l,0)/mH = 5.9754 x 10^^p{l,0). 

The disk temperature in models 1—8 is fixed to two thirds of the effective temperature 
of the central star (Hummel & Vrancken 2000), i.e., 16000 K, while different mass loss rates 
are assigned from 1.0 x 10^^^ to 2.0 x IO^^^Mq yr"^ for models 1—8 in this order. Model 
9 has the same mass loss rate as model 7, but with the disk temperature of 12000 K. This 
case is introduced for the examination of the energy balance. For models 1—8, we obtain 
the rotational velocities of v^{zu = 1) = 590 km s~^ and V(p{zu = 1000) = 5.35 km s~^ and the 
expansion velocities of v^{zu = 1) = 0.00472 km s~^ and v^(w = 1000) = 30.6 km s~^. For model 
9, we obtain V(i,{zu = 1) = 591 km s""*^, v^{w = 1000) = 6.30 km s""*^, v-^{w = 1) = 0.00303 km s~"^, 
and fro(ti7 = 1000) = 22.4 km s~^. Hence, the rotational velocities at zu = 1000 are about one 
third of the Keplerian velocity, while the Mach numbers at to = 1000 are about two. In figure 
1, we show the iso-contour of the number density on the meridional plane for our model 7, 



^ In conversion from the column density to the number density in Iwamatsu & Hirata (2007), there was an 
error that the column density is divided by /i = 0.5. Therefore, the correct density is a half of the tabulated 
one, and all of the mass loss rate appeared in Iwamatsu & Hirata (2007) should be corrected by multiplying 
two. 



Fig. 1. Iso-contour of the number density on the meridional plane in model 7. The outer boundary is 
also shown by a solid curve. 

together with the disk boundary which will be described in section 3.5.2. 



3. Non-LTE Calculation 

3.1. Statistical Equilibrium Equation 

We assume that the disk consists of only hydrogen atoms. When we designate the 
uppermost bound level of hydrogen atom by uq, the charge conservation is expressed as 

no 

N, = N+ = N{w,z)-Y,Nn , (4) 

n=l 

where Ne,N+, and A'^ are the electron number density, the proton number density, and the 
population of the n-th level at the point {w,z), respectively. The 2s and 2p sublevels are treated 
separately as in Marlborough (1969). We assume complete redistribution for the line radiation. 
Following the notation of Kriz (1974), the statistical equilibrium equations are expressed as 

no 

+ E {^n'An'n + Nn'Bn'nJnn' — NnBnn'Jnn') 
n'=n+l 

n-1 

— E i^n^nn' + NnBnn'Jn'n — ^n'Bn'nJn'n) + N^N^Un 
n'=l 

-47ry a„(z/)-^|l-— exp(^-— j|c/z/ = (n = 1,2, ■ ■ ■ ,no) , (5) 

where Juc-, Jnn' are the mean intensity of the continuum and line radiation, respectively. An'n-, 
Bnn'i and -B^'n {f^' > n) are Einstein coefficients for the spontaneous emission, absorption, 
and stimulated emission, respectively. Cnn' is the collision coefficient, and a„(z/), a„ are the 
ionization and recombination coefficients, respectively. Here, we neglect two-photon emission 
between Is and 2s levels. Introducing a departure factor 6„, = Nn/N*, where A'^ is the LTE 
value, equations (5) are written as 

71—1 n— 1 

E bn'gn'e""^' (NeCn'n) " &n [ E dn'e"""' (iVeC„'„) 

n'=l 7i.'=l 

no n-1 ^nr, J 

•Jv 



+ ^„e^" A^e( E C'nn' + C„e)+ E / a„ (i/) -^^rfz/ 

n'=n+l n'=l •^'^^ 

+ E ^n'{^ne^"(iVeC„„0 +^n'e^"'v4„,„{W}} = -2( — ^ — f 

n'=n+l 

-g^e^"{N,C^, + A7rJ^ a„(z/)-^exp ( - — )dz/} (n = 1,2, ■ ■ ■ ,no) , (6) 
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where {nn'} is the net radiative bracket for H„^„/ hne defined in Thomas & Athay (1961) and 
is given by 

{nn'} = l-^ rp^ {n<n'). (7) 

3.2. Single-Flight Escape Probability 

We approximate the net radiative bracket by a single-fiight escape probabihty (e.g., 
Kastner 1994) for easier treatment of the hne radiative process. This approximation is intro- 
duced as an extension of Sobolev approximation for the disk in which the velocity gradient is 
smaU. The single-fiight escape probability Pnn' is defined as 

I3nn' = -ri cii^?/'(z/- z/„„/)exp{-r^„,(i^)} , (8) 

where Q is the solid angle, z/„„/ is the central frequency of H„„/ line, and r^n'i^) is the optical 
thickness from a point in the disk to the disk boundary for H„„/ line at the frequency u along 
the path, for which we must take into account the velocity field. The function ip is an intrinsic 
broadening function and we assume thermal Doppler broadening. Then, ip is given by 

where Aud is the thermal Doppler width, defined as 



Aun = ^ J^^. (10) 

C V "^H 

3.3. Continuous Radiation 

The continuum mean intensity J^c is the sum of the mean intensity of the direct and 
diffuse radiation: 

Ji^c = J,^c "I" "^uc- (1-'-) 

The direct radiation from the star, J*^, is given by 

J:, = ^/ Ke-^^dQ, (12) 
Air Jn^ 

where /*, r^, are the stellar radiation intensity at the stellar surface, the continuum optical 
depth at the frequency u, and the solid angle subtended by the star, respectively. When ^ 0, 
J*^ is approximated by 

JL^wi:, (13) 

where is a geometrical dilution factor expressed by 



Here R is the distance from the center of the star in units of -R*. 
The diffuse radiation, J^^, is expressed as 



(15) 



where I^^i^) is the diffuse radiation intensity towards the direction s, which is given by the 
formal integral of the radiative transfer equation as 



Here, 3,^ is the source function and r^^ is the optical thickness from a point in the disk to the 
disk boundary. Although electron scattering yields the observed polarization, we do not include 
this effect in our calculations because it does not affect the energy conversion in the disk and 
because the scattering process makes the convergence in the iteration slower. 

3.4- Physical Constants 

a„(i^) and a„ are adopted from Burgess (1964) for n = 2s and 2p. We adopt Mihalas 
(1978) and Johnson (1972) for the other a„(z/) and a„, respectively. The gaunt factors for 
the bound-free and free-free processes are taken from Kurucz (1970). The free-free absorption 
coefficient kgi^u) and the continuous emission coefficient ji^ are adopted from Mihalas (1978) and 
Brussaard & van de Hulst (1962), respectively. An'n for n = 2s or 2p are from Karzas & Latter 
(1961), while the other v4„/,„ values are from Wiese, Smith & Glennon (1966). C„„/ and Cnc are 
obtained from equations (36) and (39) in Johnson (1972) for n,n' ^ 2s, 2p. C2s,2p is computed 
from equations (3) and (55) in Seaton (1955). Ci,2s is obtained through the interpolation of 
table X in Callaway (1985). Because 6*25,3 and C2p,z are the sum of C2s,3s, 6*25,3^, C2s,3d, and of 
C2p,3s, C*2p,3p, C2p,3d, respectively, these six coefficients are computed through the interpolation 
of table HI in Callaway et al. (1987). The coefficients C2s,n' and C2p,n' ^01 n' >4 are estimated 
from the collision rate of Johnson (1972), in the same manner as in Krolik & McKee (1978). 
For C2s,c and C2p,c, the assumption by Drake & Ulrich (1980) is adopted, i.e., Cni,c = Cnc- 

3.5. Computation 

3.5.1. Computational Scheme 

By assuming trial values for 6„, the coefficients for 6„ in equations (6) are obtained 
through integrations for the escape probability defined by equation (8) and for the mean inten- 
sities of continuous radiation. Therefore, the calculations of equations (8), (12), and (15) and 
the solutions of equations (6) and (4) constitute an iteration loop in this order. We solve the 
transfer equation for the continuous radiation by this simple integral form, though it is known 
that the convergency is very slow when the disk becomes optically thick. In actual calculations, 
we assume the uppermost bound level, uq = 10, while we adopt no = 5 when the energy balance 




(16) 
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is discussed in section 4.2. For model 1, the initial values were set to 6„, = 1/W. Then, the 
results of model 1 were adopted as the initial values for model 2, and so on. 

3.5.2. Grid points and Interpolation 

The vertical boundary of the disk is set to the surface with N{w,z) = lO^cm"^, while 
the disk radius, w^isk, is cut off at 95% of the radius with N{w,0) = 10^ cm"'^. The grid points 
for calculating the 6„ factors are set denser near the star and the equator. The fourteen grids 
are set exponentially for w, while nine grids are set parabolically for z (2; > 0), including the 
boundary. Thus, the total number of the grids is 126. We set the grid points in the following 
form 

^* = ^disk = 1, 2, ■ ■ ■ , 14), and Zij = Zi^s (^0 (j = 0, 1, ■ ■ ■ , 8) , (17) 

where Zi^s is the vertical boundary with N{wi,z) = lO^cm"^. Hereafter, we call these grid points 
as population grid points. Analytic setting of the grids is convenient for the interpolation in the 
{i,j) space with an equal interval. Interpolation is made only for the bn factors, since all the 
relevant quantities are calculable through the 6„ factors in the isothermal disk. Actually, the 
log6„ values are interpolated. Easy and stable interpolation is required for precise integrations 
described below. 

3. 5. 3. Integration 

Equations (6) include many integrations in the coefficients of 6„. It is evident that 
the iteration converges to the false values if the integrations appeared in equations (6) are 
inaccurate. The 1% accuracy of the integration is aimed in our iterative calculations and also 
in calculations of the emergent flux and the line profile. We checked the accuracy for several test 
integrands in each integrations appeared in equations (6). Typical test integrands are based on 
the assumption of 6„ = 1/W. The integration is made mostly by Simpson's 3/8 formula. 

The optical depth r^„/(z^) in equation (8) is calculated with the meshpoints concen- 
trated around the place where the relative velocity is zero along the path. The meshpoints 
for frequency integrals are set with an interval of 0.08Az/£) for the range of — 5Az/£) to 5Az//). 
The angular integration in equation (8) is performed on a polar coordinate with a polar axis 
parallel to the rotational axis of the star. The meridional integraion is performed with 25 mesh- 
points, which distribute denser near the equator, and 24 meshpoints are set for the azimuthal 
integration with an equal interval. 

The optical depth in equation (12) is calculated with 120 meshpoints concentrated 
in the regions where the function N'^{w^z)/W is large. The angular integration is performed 
by introducing a polar coordinate whose polar axis passes the center of the star. We adopt 18 
meshpoints from to 9^ (the angular radius of the star) for the meridional integration, and 36 
meshpoints from to 27r for the azimuthal integration. 

Since the continuous emission produced in the denser part near the star is effective for 
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the diffuse radiation, J^^, in equation (15), we introduce a polar coordinate whose polar axis 
passes the footpoint, {vj,z) = (1,0) for the angular integration. The angular meshpoints are 
set denser near the footpoint and towards the equatorial plane, depending on the population 
grid point. The meridional integration is performed first, and then the azimuthal one. In the 
case of the population grid point distant from the star, about 100 meshpoints are required for 
the meridional integrals, while 12 meshpoints is enough to satisfy the 1% precision when the 
population grid point is located in the dense part, due to its isotopic nature. For the azimuthal 
integration, 26 meshpoints are set, which are distributed denser near the equator. 

The path integral in equation (16) is made by using Simpson's 3/8 rule for the calculation 
of the optical depth, while the trapezoidal rule is applied for the integration of the source 
function. The segments for the calculation of the optical depth are divided repeatedly until they 
have satisfied the constraint that the absorption coefficient ratio for both edges of the segment 
is between 1/1.2 and 1.2 for all frequency meshpoints. Moreover, the additional constraint is 
added for the optical depth of the segment, At^,, at frequency u: the segment is divided until 
satisfying At^ < 0.2 for < < 2, Ar^ < 2 for 2 < < 5, Ar^ < 5 for 5 < < 10, and we 
neglect the contribution from the more distant region with > 10. As a result, the number of 
the meshpoints along the path which passes only the optically thin region is several tens, while 
it reaches several thousands when the path includes the optically thick region. 

Totally sixty four frequency meshpoints are set for the integration over the frequency 
to calculate the ionization rate in equation (6). They consist of 13 meshpoints for the Lyman 
continuum, 10 meshpoints for the Balmer continuum, 7 meshpoints for the Pashen, Brackett, 
and Pfund continua, and 4 points for the other continua with n > 5. 

3.5.4- Convergency 

Our simple iteration scheme requires many iterations in the case of optically thick disk. 
In such a case, the convergence is slowest for bi at the place slightly distant from the star 
around the equatorial plane where the Lyman optical thickness towards the disk boundary is 
very large. When the ratio of the new bi value to its previous one keeps monotonous variation 
during five consecutive iterations, the five power of the ratio is multiplied to the new bi factor 
in order to accelerate the convergence, while the root of previous and new 6„ factors is adopted 
as a next trial value in order to relax the solution. 

By examining the convergence behavior of 6„, we found that all the 6„ factors except 
bi in the dense region converge when all the ratios of newly calculated departure factors 6i 
to the previous factors are within 1 ±0.01 at all the population grid points. Such a behavior 
in iteration is understandable since our models are not optically thick towards the vertical 
direction in the Balmer and the other continua with longer wavelength, and since the fully 
ionized state is assured from the assumed high disk temperature. As will be shown in section 
4.3, Lyman continuum does not play so important role in the radiation field in the disk. Then, 
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we decided to stop the iteration when the above 1% criterion is fulfilled for all 6„ values. When 
we compare such a result with the results from further iteration and from the iteration with 
the initial bi values in opposite sense, a few tens percents deviation at most from the guessed 
convergent values was noticed for bi in the dense region of model 7. 

3.5.5. Computational Time 

All calculations in this paper were made with the computer with the CPU of the AMD 
Opteron 250 (2.4 GHz) at the Astronomical Data Analysis Center, ADAC, of the National 
Astronomical Observatory of Japan. The calculation times for one iteration are about twelve 
hours for the optically thin disk, about eighteen hours for the optically thick disk. The 1% 
criterion was achieved after ten or less iterations for the optically thin case, and after about 
forty iterations for the optically thick case. The total time is approximately several days for 
the optically thin case and two weeks for the optically thick case. 

4. Global Features 

4.1. Non-LTE State 

In order to show the optical property of our model sequence, we list, in table 2, the 
respective continuum optical thicknesses, Tr,Tz, and at each continuum limit from the foot- 
point [w = l,z = 0) to the disk boundary towards the radial, vertical, and azimuthal directions 
for our models 1, 4, and 7. Reminding an observational fact that the second Balmer jump due 
to the strong Balmer shell lines appears in strong shell stars (Divan 1979), we cover the range 
of the mass loss rate for which the radial optical thickness at the Balmer shortward edge is 
larger than unity. Model 1 is optically thin, while model 7 is optically thick at all shortward 
edges. In the wavelength where the optical depth for the continuum is large, the role of diffuse 
radiation becomes important in the non-LTE state of the disk. In model 1, the contributions 
of the direct and diffuse radiation to the ionization from the ground level are comparable in 
the vicinity of the star, while the direct stellar radiation overwhelms in the longer wavelength 
region. In model 7, the diffuse radiation overwhelms the direct radiation in the dense part of 
the disk. Taking into consideration that our BIV star emits the radiation in the Lyman and 
Balmer continua by 0.2% and 93.3% of the total luminosity, respectively, we classify the optical 
property of the disk by the optical thickness for the Balmer continuum. We regard that models 
1 and 7 are the respective typical cases for the optically thin and thick disks. Notice that the 
optical thickness at the shortward edges of Balmer and Paschen series has the smallest values 
among those at all the shortward edges. These facts imply that the stellar Balmer continuum 
principally governs the radiation field in Be stars. We also note that an inequality < < 
holds not only at the footpoint but also in the main part along the equator. 

Figure 2 shows the variation of the bn factors along the equatorial plane as a function 
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Table 2. Continuum optical thickness towards the radial direction (rr), the vertical direction (tz) and the azimuthal 
direction (r^) from the footpoint to the boundary in our three models. All the bound-free processes and the free-free 
process are included. The first row corresponds to the threshold wavelength: the integer means the bound level 
related to the bound-free transition, while the ' — ' and signs indicate the shortward and the longward of the 
threshold, respectively. For example, 'Ic— ' means the shortward edge of the Lyman limit. a{b) means a X 10*. 



model 1 model 4 model 7 
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of the distance for our models 1, 4, and 7. The geometric dilution factor W is also shown with 
the dashed lines. The 6„ factors once increase outwards from 6„ ~ 1 in all our models, because 
the available stellar radiation decreases outwards. The overpopulation is less remarkable in the 
higher levels, since the coUisional process operates more effectively in the higher levels. The Is 
and 2s departure factors increase towards the outer boundary, because they are the ground level 
or the metastable level. It is noteworthy that the second level population approaches to 1/W, 
but does not exceed 1/W, reflecting a fact that the stellar Balmer continuum is a principal 
source of ionization. While, the 6„ factors of the excited levels (n = 2p, 3,4. . .) turn to decrease 
at some distance, and finally converge to the nebular case A. Similar behavior was concluded by 
Poeckert & Marlborough (1978) (see figure 2 of Hirata & Kogure 1984). We confirmed that the 
turning point with maximum 6„ reflects the change from optically thick nature to thin nature in 
H„_i „ line, i.e., the leading member of the series, thus degrading the n-th level population. The 
turning point shifts outwards in the denser disk with larger mass loss rate, caused by optically 
thicker nature. In the outer region, the disk becomes optically thin for the lines, first, from 
the higher series members, then, to the lower series members. In the optically thick disk, the 
LTE region is added in the innermost part where it is optically thick in all wavelength regions. 
Figure 3 shows the iso-contour map of the departure factors 6.„ on the meridional plane for our 
standard model 7 with uq = 5. The 6„ factors become smaller as leaving from the equator, 
reflecting the smaller optical thickness and subsequent more important role of the direct stellar 
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Fig. 2. bn factor variation on the equatorial plane for the models 1, 4, and 7. The dashed curves indicate 
1/W. 



Fig. 3. Contour map of the departure factor 6„, where n = l,2s,2p,3,4,5 on the meridional plane for model 
7. The X signs indicate the position of the population grid points in our computation. The iso-contours 
of the dilution factor W are also shown in each panel. 

radiation from the polar region. For comparison, we also show the contours for some 

typical values. 

4-2. Energy Balance 

Next, we calculate the radiation energy gain and loss in our isothermal disk and examine 
the deviation from the radiative equilibrium. The calculations were made for rio = 5. The energy 
loss at any point in the disk is given by 



n'=2 n=l 

where the first and second terms correspond to the loss through the continuous process and the 
line process, respectively. The energy gain is expressed as 



If the energy loss is larger than the energy gain, the disk temperature should become 
lower than the assumed one when the radiative equilibrium holds. We define the ratio of 
^^loss to -Egain by Rig. Figure 4 shows the ratio, -Rig, along the equatorial plane for models 1, 
7, and 9. The model 9 (Tjisk = 12000 K) is shown for illustrating the case of the lower disk 
temperature. The label T represents -Rig, while the labels C and L indicate the continuum and 
line contributions, respectively. 

In the case of optically thin disk (model 1), our isothermal assumption with Tdisk = 2/3Tcfr 
fills almost the radiative equilibrium on the equatorial plane. In model 7, however, the energy 
loss is larger than the energy gain in almost all part of the equator. There exist two regions 
where the value of -Rig reaches maximum: the inner part around 4R^, and the outer part around 
60-R*. The formation of the inner cool region is caused by the ultraviolet radiation deficiency 
due to the optically thick nature. The outer cool region is formed because Lyman a becomes 
optically thin, thus the La cooling becomes effective. Outside of the La cooling region, the 
factor drops drastically, so the La cooling becomes ineffective. The energy loss by the line 
process is negligibly small in the inner part, while it reaches a half of the total loss in the outer 
La cooling region. This La cooling region will not be observed in the optical spectral region, 
because no contribution is expected due to its low density. In model 9 (Tdisk = 12000 K), the 
global feature resembles to that in model 7. Note that the disk temperature in the innermost 
LTE region should be higher than the assumed one, in the sense that it approaches to the 




no n' —1 



(18) 




(19) 
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Fig. 4. Ratio of tiie energy loss to tiie energy gain along the equatorial plane. The label T represents the 
ratio for the total energy loss R\g, while the labels C and L indicate the continuum and line contributions, 
respectively. 

Fig. 5. Ratio of the energy gain rate to the loss rate on the meridional plane in models 1, 7, and 9. The 
upper three panels display the region within w < 20, while the lower panels include more distant region 
up to < 170. 

stellar effective temperature. 

Figure 5 shows the contour maps of -Rig on the meridional plane for models 1, 7, and 9. 
One can see that the cooler region is concentrated to the equatorial plane and the disk is hotter 
in general as 1^1 becomes larger, caused by the penetration of the stellar ultraviolet radiation. 

4-3. Energy Flow and Conversion 

Let us discuss the radiation flow and conversion in the disk. Since the continuum op- 
tical thickness is smallest towards the vertical direction (see table 2), the continuous radia- 
tion emitted in the disk tends to escape from the upper and lower surfaces of the disk. The 
inclination- angle dependency of the emergent flux will be given in the next section. 

Although the La cooling is effective in the outer part of the disk, it is found that the 
total amount of energy emitted in the lines from the whole disk is negligible, compared with 
that in the continuum. This is because the amount of the energy emitted in the region far 
from the star is much smaller than that near the star. Hence, we consider only the continuous 
radiation in the following. We define the continuum luminosity in Lyman, Balmer, ■ ■ ■ continua 
(Li, L2, ■ ■ ■ ;n = 1 — 6) for rio = 5 by 



where Fi,{i) is the flux integrated over the disk and the star towards the inclination angle i at 
the frequency z/, z/„ being the threshold frequency of each continuum band with the exceptions 
of z/q = 00 and i^e = 0. We also define the similar quantity L* for the stellar radiation. 

In figure 6, we show, for models 1, 7, and 9, the difference L„ — at each band, as 
well as the difference of the total luminosity, L — L^ (designated as "all"), all being normalized 
by the stellar bolometric luminosity L^. The total luminosity of model 1 is slightly less than 
the stellar luminosity, indicating the slightly higher isothermal equilibrium temperature. Total 
luminosity of model 7 is larger than the stellar luminosity, while that of model 9 is less than the 
stellar luminosity. This indicates that the isothermal equilibrium temperature is in between 
12000 K and 16000 K for model 7. This conclusion is consistent with the results in the radiative 
equilibrium models calculated by Millar & Marlborough (1999), Carciofi & Bjorkman (2006), 
and Sigut & Jones (2007), since they all gave Tdisk around 0.6Teff as an equivalent isothermal 
disk temperature. 




(20) 
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Fig. 6. Disk luminosity in each band relative to the stellar bolometric luminosity for models 1, 7, and 9. 
The abscissa expresses each band, e.g., lc=Lyman continuum and 'all' indicates the total luminosity. The 
scale of the ordinate in the upper and lower two panels are different. 

Fig. 7. Spectral energy distribution for three models. 

Next, we discuss the radiation energy conversion in the disk. It is well known that stellar 
ultraviolet radiation is converted into radiation with longer wavelength and line radiation in 
the Be star disk. Figure 6 clearly shows that Balmer continuum photons are converted into 
photons with longer wavelength. Similar conclusion was reported by Millar & Marlborough 
(1999). Stellar Lyman continuum plays little role in the disk except the ground level regulation 
in the outer envelope of the disk where the optical depth for Lyman continuum towards the 
central star is small. In the optically thick disk, only a small fraction of the Lyman continuum 
is converted into the optical and infrared radiation (models 7 and 9). While, in the optically 
thin disk (model 1), Lyman continuum photons are even produced in the disk from the Balmer 
continuum photons. This result comes from a fact that a Bl V star emits almost all the radiation 
not in the Lyman continuum but in the Balmer continuum. This is a general conclusion for 
Be stars, because the stellar Lyman continuum is more deficient in Be stars with later spectral 
types. 

5. Continuum Property 

5.1. Spectral Energy Distribution 

We describe several observable quantities for our models with T^isk = 16000 K in present 
and next sections. Figure 7 shows the spectral energy distribution (SED) for models 1, 4, and 
7. The SEDs for different inclination angles are shown, together with the stellar SED. The 
infrared excess develops as the mass loss rate increases and as the inclination angle decreases. 
The Lyman continuum shows the deficiency or excess, depending on the inclination angle in 
the case of optically thick disk. In model 7, the disk continuum contributes to the SEDs even 
in the Balmer and Paschen continua. 

5.2. UBV Bands 

The UBV magnitudes were calculated by using the response functions for the photo- 
metric filters given in Buser & Kurucz (1978). Figure 8 shows the V magnitude variation as 
a function of the mass loss rate. The brightness in the V band increases as the mass loss 
rate increases and the inclination angle decreases. The brightness decrease in the V band is 
not resulted and only slight brightening in the V band is obtained in the equator-on case in 
our present model sequence (Tdisk = 16000 K). We note, however, that the V band brightness 
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Fig. 8. V^-band magnitude variation as a function of the mass loss rate for the model sequence with 
Tdisk — 16000 K. We also show the disk mass in the upper part. 

Fig. 9. Flux variation as a function of cosi in the UBV bands for model 7. 

decrease by 0.2 mag is obtained for i = 90° in model 9 (Tdisk = 12000 K). The observed range 
of AV" ~ 0.5 mag is well covered by the mass loss rate less than 10~^° Mq yr~^. 

Next, we examine the inclination angle dependency of the emergent flux in each band. 
Figure 9 shows the fluxes in the UBV bands for model 7 as a function of cosz, to which 
the projected area is proportional. It is seen that the emergent fluxes are, roughly speaking, 
proportional to cost, though they are slightly concentrated towards the lower inclination angle. 

Figure 10 shows the variations on the color-magnitude (CM) and color-color (CC) dia- 
grams. The B — V color becomes redder and the brightness in the V band becomes larger as 
the mass loss rate increases. The amount of variation is larger for smaller inclination angle. 
It is interesting to note that the variation on the CM diagram almost lies on the same locus, 
regardless the mass loss rate and the inclination angle. The color-color plot in figure 10b indi- 
cates that the U — B color becomes bluer for the low inclination angle and redder for the high 
inclination angle, as the disk develops. We confirmed that no color change in U — B occurs at 
i ~ 70°. We note that our present results are consistent with the observed variations on the 
CM and CC diagrams for early Be stars (e.g.,Hirata & Hubert 1981). 

5. 3. Far-Infrared Region 

The free-free transition is a dominant process in this wavelength region. Therefore, the 
emergent flux is rather free from the non-LTE state and reflects the density and temperature 
structure of the disk more directly. Figure 11 shows the magnitude variation at four IRAS 
wavelengths as a function of the mass loss rate. The monochromatic magnitudes were calculated 
at 12, 25, 60, and lOO/xm, since the SED in this wavelength region is not complex. The 
brightening in magnitude scale increases linearly with the increase of the mass loss rate in the 
pole-on case. 

The difference between z = 0° and 90° reaches 1.5—1.8 mag in the case of the larger mass 
loss rate, resulting one order of magnitude difference in the mass loss rate in our models. We 
suggest that the pole-on approximation in Waters (1986) is valid for i < 30° and the systematic 
deviation is expected for the higher inclination angle. Similar result has been concluded by 
Carciofi et al. (2007). We found such a correlation between the mass loss rate given by Waters, 
Cote & Lamers (1987) and the projected rotational velocity vsini given by Cote & Waters 
(1987), as is shown in figure 12. There is a tendency that the mass loss rate becomes larger for 

Fig. 10. (a) Variation on the Color-Magnitude diagram, (b) Variation on the Color-Color diagram. 
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Fig. 11. Monochromatic magnitude change at four IRAS wavelengths (12, 25, 60, and 100 /im). 

Fig. 12. Mass loss rate given in Waters, Cote & Lamers (1987) against usini. Different spectral types 
are distinguished by different symbols. 

larger stellar projected rotational velocity vsini. We also note that, in figure 12, the stars with 
large vsini (~400 kms~^) have the mass loss rates one order of magnitude larger than that 
with the smaller vsini (~ 100 kms~^) for all the spectral type groups. Our rough estimate 
suggests that the mass loss rate based on our model is 10~^^ — 10~^° MQyr~^, three orders of 
magnitude smaller than that of Waters, Cote & Lamers (1987), reflecting the difference in the 
initial velocity 1^^(117 = 1). 

Figure 13 shows the flux variation at four IRAS wavelengths as a function of cosi for 
model 7. Note that the disk contribution overwhelms the stellar flux in the IRAS wavelength 
region. Again, the flux increases almost linearly with cosi. However, due to the opaque nature 
in this wavelength region, the disk has a finite opaque area, thus resulting the infrared excess 
even in the equator-on view. Note that the cost dependence of the infrared flux is convex, in 
contrast to that in the optical region (see figure 9). 

6. Line Property 

6.1. B aimer Line Profiles 

Balmer line profiles were calculated under the assumption that the intrinsic broadening 
function is given by the thermal motion to see the effect of velocity field, although Stark 
broadening becomes important in the denser part of the optically thick disk. Our models 
cover the large range of the continuous optical thickness at Ha from = 0.01 (model 1) to 
Tr = 375 (model 8) (see also table 2). In the optically thick models, the contribution from 
the disk continuum is noticeable even in the visual region, as is shown in the previous section. 
This means that the resulting Balmer line profiles are complex, disturbed not only by the 
shell absorption, that is, the line absorption of the photospheric radiation, but also by the line 
absorption of the disk continuum. The latter is attributable to the so-called pseudo-photosphere 
(Harmanec 1994, 2000). 

In figure 14 we show the rectified Ha, H/3, and H7 line profiles in three models 1, 4 and 
7 for various inclination angles. The photospheric component is not included in this figure. In 
the case of i ~ 0°, the narrow double-peak profile appears. This special case will be discussed 
in section 6.2. The wine-bottle type profile is traced for i = 15° and 30° in model 7. The broad 
double-peak profile develops for the larger inclination angle. The peak separation is larger 

Fig. 13. Flux variation as a function of cosi in the IRAS bands for model 7. 
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Fig. 14. Ha, H/3, and H7 line profiles for our three models. The photospheric line is not considered in 
these profiles. The dotted, dashed and solid lines are the line profiles for models 1, 4 and 7, respectively. 
Note that the top panels have the vertical scales different from the other lower panels. 

Fig. 15. Ha line profiles with the inclusion of the photospheric component in model 4. The contribution 
of the photospheric line, the shell absorption and the disk emission component are also shown. The 
ordinate is in the flux form. 

for the smaller mass loss rate and for the larger inclination angle. The formation of the peak 
separation will be also discussed in section 6.2. The strong shell component appears for i = 90° 
in models 4 and 7. These behaviors, including the narrow double-peak emission line profile 
in the pole-on case, are similar to those obtained by Hummel (1994). We also note that the 
central quasi-emission (CQE) is noticed for i = 90° of model 1 (for CQE, see Hanuschik 1995). 

The Ha line profiles, in the flux form, with the inclusion of the photospheric compo- 
nent are illustrated in figure 15 for model 4, where we also show the contribution from various 
components: the emission component, the shell absorption component, and the photospheric 
component. The emission component includes the line and continuum originated in the disk. 
The intrinsic photospheric profile was taken from ATLAS9 (Kurucz 1993) with the same pa- 
rameters as for the stellar SED. This profile is put on the rotating stellar surface with the 
critical equatorial rotational velocity of 591 kms~^ without limb darkening. One can see the 
veiling effect by the disk continuum, which is more noticeable in the smaller inclination angle. 

For illustrating a case contaminated heavily by the disk continuum, we show, in figure 
16, the Ha profiles in model 7. In addition to the veiling effect, the emission component has 
absorption features below the disk continuum in the wing part for i = 30° — 75° and in the 
central part for i = 90°. These absorption features are formed by the line absorption of the 
disk continuum, i.e., in the pseudo-photosphere. This suggests that 1) for the intermediate 
inclination angle, it is difficult to segregate the photospheric component from the observed 
profile, caused by the broad shell absorption and the pseudo-photospheric absorption, and 2) 
for i ~ 90°, the observed central intensity could be almost zero even if the disk continuum 
is present. We note that the emission equivalent width rectified by the disk continuum and 
integrated over the line feature may become even negative, caused by the dominant pseudo- 
photospheric absorption feature. This is notable in the higher Balmer lines in our models. We 
also note that a weak central shell absorption feature appears for i > 75° in model 7. Our 
major results on Ha line are similar to those of Hummel (1994) except the veiling effect and 
the pseudo-photospheric absorption, which are newly introduced in this paper. 

Fig. 16. Same as figure 15 but in model 7. 
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Fig. 17. Vertical Ha source function variation at w — 1.06 in model 7. Several optical depths at the line 
center are indicated. 

6.2. Peak Separation 

As stated in section 6.1, the narrow double-peak emission profile appears in the case of 
2 ~ 0°. The same characteristics was concluded and interpreted by Hummel (1994) in terms of 
the non-coherent scattering (NCS). Figure 17 illustrates the vertical change of the Ha source 
function at the distance of = 1.06 in model 7. The optical depth at the line center is also 
shown. One can see that the source function increases as the optical depth increases in the 
effective line formation region. Similar tendency is kept in the inner part of the disk, where the 
source function is large. Recalling that no velocity gradient exists in the vertical direction in 
our model, the central reversal is naturally formed. 

In figure 18a, we show the resulting peak separation as a function of sini. As sini 
increases, the profile characteristic changes from double-peak profile caused by NCS, passing 
the single-peak profile, to the double-peak profile of rotation origin. The peak separation caused 
by NCS is within two or three Doppler widths in our calculations and increases as the mass loss 
rate increases, caused by the stronger central reversal. The domain of this type of double-peak 
profiles expands further to the larger inclination angle as the mass loss rate increases. The 
same type profiles are also formed in H/3 and H7 lines with the peak separations narrower than 
those in Ha line (see figure 14), reflecting their smaller optical thickness. 

The peak separation originated in the rotational velocity field may be broadened by the 
shell and the other absorption components. In figure 18a, we also show the peak separations 
derived from the emission component only, the profile including the shell absorption component, 
and the profile including further the photospheric component for models 4 an 7. While, we 
show only the peak separation derived from the emission component for model 1, because the 
the emission component is very weak, so the determination of peak position is difficult when the 
shell absorption and the photospheric component are taken into account. One can see that the 
effect of making peak separation larger by such absorption features is small even in the weak 
emission profiles for i ~ 90° in model 4. The peak separations in model 7 is little influenced by 
the background absorption features, since the pure emission profiles have sharper peaks. Note 
that the large deviation of the peak separation at z = 90° from those at the other inclination 
angles is originated in the pure emission component itself. 

Now, let us examine the peak separation of the double-peak profile in the higher incli- 
nation angle. The peak separation of the double-peak emission line profile in Be stars has been 
often used for estimating the outer radius of the disk. Observationally, it is well known that 
the peak separation is smallest in Ha and becomes larger for the higher Balmer lines, reflecting 
the emitting sizes in each lines. The outer radius of the disk, -Rdisk, in the Keplerian rotation 
law is estimated as 
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where AV^eak is the peak separation of the double peak profile. This formula is derived from 
the profile calculation for the rotating flat disk with a finite outer radius, under the assumption 
that the line surface brightness distribution is given (e.g., Huang 1972, 1973, Hirata & Kogure 
1984). The velocity at the maximum intensity corresponds to that at which the iso-line-of- 
sight velocity curve contacts with the outer radius, unless the surface brightness distribution is 
a steeply decreasing function of the distance. If the surface brightness is locally proportional 
to tu""*, it is easy to show that the region with m > 2.5 does little contribute to the emission 
profile in the Keplerian disk (see equations (4) — (7) in Hirata & Kogure 1984). Then, the peak 
separation in our model which extends infinitely towards the radial direction is thought to 
reflect the maximum radius with m < 2.5. 

The peak separation given by equation (21) is proportional to sinz when the outer radius 
is fixed. While, in figure 18a, we notice a quasi-linearity except for i ~ 0° and i ~ 90°. Hence, 
the calculated peak separation reflects somehow the effective disk radius, although the quasi- 
linear part does not pass the zero point, in contrast to the implication from equation (21). It is 
possible to derive the effective outer radius by fitting the slope in the quasi-linear part of figure 
18a to that in equation (21). Then, we obtain -Rdisk = 1-7, 4.9 and 12.3 for models 1, 4, and 
7, respectively. 

Although we did not calculate the line surface brightness for our models, it is possible 
to examine how the peak position is formed in our models. Figure 19 shows the variation of 
the escape probability (3, the source function S, and p = Ai'3X]n=2s,2p^(3,'T-)/5(?^,3) for Ha line 
in model 7 along the equatorial plane. The last quantity is proportional to the Ha photon 
energy which escapes from the equatorial plane to the outside of the disk and may be regarded 
as the angle-averaged surface brightness distribution on the equatorial plane. With a slope 
of m = 2.5 for p, we obtain -Rdisk = 1-3, 3.5, and 7.8 -R* for models 1, 4, and 7, respectively. 
These values are, roughly speaking, consistent with those estimated from the slopes in figure 
18a, mentioned before. The source function illustrated in figure 19 essentially behaves in the 
same manner as p in the main part of the disk, having a slower decline in the inner part and a 
steeper decline in the outer part. Hence, we can guess the formation region of any lines from 
their source function behavior along the equatorial plane. The source function for H„„/ line is 
approximately proportional to 6„/ /6„. From figure 2, it is evident that the higher Balmer line is 
formed in the region nearer to the star than the region where the lower Balmer line is formed, 
thus has a larger peak separation. 

In figure 18b, we show the Ha outer radii estimated from equation (21), together with 
the radii estimated fromp and /5(2p,3) = 0.5 at the equatorial plane. We choose P{2p,3) = 0.5 as 
the effective boundary between the optically thick and thin regions at Ha. Roughly speaking, 
the radii estimated by equation (21) fall in between our two rough estimates. However, the 
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Fig. 18. (a) Ha peak separation against sini for tliree models. Tlie peak separation obtained from the 
emission component, the emission+shell component, and the emission+shell+photospheric components 
are designated by "+", o, and x, respectively, (b) Disk radius calculated from equation (21). Symbols 
are the same as in (a). We also show two rough theoretical estimations: the distance where the escaping 
Ha photons begins to decrease rapidly from the equatorial plane (thin dotted line with a designation of 
m = 2.5) and the distance where the escape probability of the Ha photon becomes 0.5 at the equatorial 
plane (thick dashed line with a designation of /3 = 0.5). Sec the text for more details. 

Fig. 19. Escape probability, source hmction, and p = ^3j2n=2s 2p ^('^i ")'^("'' 3) equator for Ha 

line in model 7. Sec the text for the details. 

direct application of equation (21) to the observation gives only a rough estimation of the 
disk radius. Hummel (1994) obtained the same conclusion and discussed the peak separation 
formation, based on the emissivity behavior for his A^(ci7,0) oc models. 

6. 3. Equivalent Width and Balmer Decrement 

Physically meaningful line emission strength would be the total flux calculated from the 
emission component above the disk continuum, illustrated in figures 15 and 16. We call it a line 
emission power, hereafter. However, the synthesized emission profiles are more complex one, 
disturbed by the disk continuum, by the shell component in the case of large inclination angle, 
and by the pseudo-photospheric absorption in the optically thick disk. Here, we simply define 
the emission equivalent width, Wq,, as the area above the apparent continuum level in figure 
14. This corresponds approximately to the equivalent width above the photospheric absorption 
profile in actual observation. 

Figure 20 shows the variation of Wa as a function of the mass loss rate. The computed 
Ha equivalent widths are rather small and are similar to those reported by Hummel (1994) 
and Hummel & Hanuschik (1997) for their models with N{w,Q) oc w~'^. They are 6 A at most 
and do not reproduce the observed ones (e.g.,Briot 1971; Dachs, Hanuschik & Kaiser 1986; 
Dachs, Rohe & Loose 1990, 1992). Sigut & Jones (2007) have also concluded that the Ha 
luminosity for their models with p oc w~^'^ is too faint when compared with the observations. 
This problem will be discussed in the last section. At the high end of the mass loss rate larger 
than 5 x 10~^^AfQyr^\ the equivalent width gradually decreases due to the veiling effect for 
i < 60°. In figure 20, we also show the equivalent width integrating over the profile, including the 
absorption feature in the case oi i = 90° (dotted line) in order to demonstrate the effect of the 
absorption feature: the broad shell absorption component overwhelms the emission component 
when the mass loss rate is low (see figures 14 and 15). This implies that the development of a 
broad absorption feature is the first spectroscopic evidence of a new mass loss event in the case 
of i ~ 90°. 

Figure 21 shows (a) the line emission power and (b) the equivalent width of Ha against 
cosi The Ha emission power is free from the veiling effect and from the shell absorption by its 
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Fig. 20. Variation of the Ha emission equivalent width as a function of the mass loss rate. The disk mass 
is also shown in the upper part. The dotted line corresponds to the equivalent width with the inclusion 
of the absorption feature below the apparent continuum level. 

Fig. 21. (a) Emission power, (b) emission equivalent width of Ha as a function of cos?. 

definition. Declines at i = 90° in models 4 and 7 are caused by the central pseudo-photospheric 
absorption in figure 21a and by the central shell and pseudo-photospheric absorption in figure 
21b. The Ha emission equivalent width in figure 21b suffers strong veiling effect for i < 60° 
in model 7. From figure 21a, we conclude that Ha photons also tend to escape towards the 
polar direction, in spite of the existence of velocity gradient. However, such a tendency is less 
remarkable, compared with the continuous radiation (see figures 9, 13, and 21a). We confirmed 
that H/3 and H7 lines show the same tendency. 

The Balmer decrement may be defined as the ratio of the fiuxes above the apparent 
continuum level in figures 15 and 16, corresponding to our definition of the emission equivalent 
width. Here, we call it observable Balmer decrement. We also define theoretical Balmer decre- 
ment derived from the line emission powers. Figure 22 shows (a) theoretical and (b) observable 
Balmer decrements in the form of the ratio of Ha to H/3 (-D34) versus the ratio of H7 to H/3 
(D54). For the observable Balmer decrement, we do not adopt the cases in which is less 
than 0.1 A, because the results are less accurate and are also meaningless from the observational 
point of view. Several cases for high inclination angles are lost in figure 22b due to this cut-off 
in H7 and H/3 lines. For comparison, we also show a nebular case B solution for Te = 10000 K 
(Brocklehurst 1971) by a large plus marked "Case B". The theoretical Balmer decrements are 
located in the left of case B, while the observable Balmer decrements extend towards right and 
lower direction, caused by the shell absorption effect for high inclination angle (see the cases 
for i = 75°). Our results in figure 22b fall in the observed area given in Dachs, Rohe & Loose 
(1990). In both diagrams, the Balmer decrement is steeper for higher inclination angle. This 
is consistent with the observational implication of steeper decrement for shell stars suggested 
by Dachs, Rohe & Loose (1990). 

7. Summary and Concluding Remarks 

The non-LTE state of hydrogen atoms in the isothermal steady-state decretion disk 
around a BIV star has been calculated, based on the A iteration for the continuous radiation 
and the single-fiight escape probability for the line radiation. Thus, we obtained the continuous 

Fig. 22. Balmer decrement for Ha— H7 lines. D34 is the ratio of Ha to H/? and D54 is the ratio of H7 
to H/3. (a) theoretical diagram derived from the line emission power, (b) observable diagram derived from 
the equivalent width above the apparent continuum. See the text for the details. The nebular case B is 
also indicated by a large plus sign. 
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radiation field without approximation. The basic physical process in our models was presented. 
We described the non-LTE state in terms of the bn factors. We examined the radiation flow and 
radiation conversion in the disk. The local and global deviations from the radiative equilibrium 
in our iso-thermal models were examined through counting the energy gain and loss. We also 
calculated several photometric and spectroscopic quantities and examined their dependency on 
the mass loss rate and the inclination angle. Some implications to the observable quantities 
were also presented. 

Our main conclusions are summarized as follows: 

1. The bn factors increase outwards from 6„ ~ 1. The 6„ factors of 2p and higher levels then 
turn to decrease at the distance where the H„_„_i line becomes optically thin and finally 
tend to a nebular case A solution. In the optically thick disk, the LTE region is added in 
the innermost region of the disk. 

2. Stellar Balmer radiation plays a definitive role in the radiation field of the disk. The 
conversion of the stellar Balmer radiation into the longer wavelength continuous radiation 
and the line emission is the main process in radiation field of the Be star disk. Stellar 
Lyman continuum contributes only to the ionization in the outer envelope of the disk 
where the disk is optically thin towards the central star. 

3. The radiative equilibrium almost holds in the disk with T^isk = 2/3Tefj in the optically thin 
case for the Balmer continuum. Similar conclusion is obtained by Carciofi & Bjorkman 
(2006). When the disk becomes optically thick, there appear two cool regions: the inner 
region caused by the deficient Balmer continuous radiation and the outer region caused 
by the Lyman a cooling. 

4. The optical thickness of the disk is smallest towards the vertical direction for the con- 
tinuum. Hence, the continuous radiation tends to escape from the disk vertically. Our 
models with the mass loss rate less than 10~^^ MqJt^^ cover well the observed variations 
in the UBV and IRAS bands. The infrared excess in the IRAS band for 90° is smaller 
than that for i = 0° by L5— 1.8 mag. We suggest that the pole-on approximation in 
Waters (1986) is valid only for i < 30°. The infrared mass loss rate based on our model 
is 10~^^ — lO"^'' M0yr~^, three order of magnitude smaller than that of Waters, Cote & 
Lamers (1987). 

5. General characteristics of the observed profiles of the classical Be stars like the wine-bottle, 
double peak, and shell profiles are reproduced in our non-LTE model, just as was shown 
by Hummel (1994). However, the isothermal steady-state decretion disk does not give the 
observed large Ha emission intensity. 

6. The pseudo-photospheric absorption feature through the line absorption of the disk con- 
tinuum appears in the case of the large mass loss rate. 

7. The peak separation of the double-peak profile is almost determined from the emission 
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profile, that is, the shell absorption and the photospheric absorption affect the peak po- 
sition only a little bit. Though the peak separation gives a rough measure of the outer 
radius of the disk, the resultant radius still depends on the inclination angle. This is just 
the same conclusion given by Hummel (1994). 

Most serious result in our calculation based on the N{u7,0) oc w~'^'^ law is that we did 
not get the large Ha emission equivalent width observed in early Be stars. As the mass loss rate 
increases, our isothermal model sequence with Tjisk = 2/3Tefj deviates more from the radiative 
equilibrium in the sense that the radiation loss from the disk is larger than the energy input 
from the star. This means that the radiation energy lost from the disk comes partially from 
the compensation of the internal energy in the disk with the temperature higher than the 
equilibrium temperature. Hence, the Ha intensity will not increase even if we introduce the 
radiative equilibrium model with the same density structure. 

It is highly probable that we should introduce some factors which enhance the line emis- 
sion, compared with the continuous emission. The Be star disk may extend more towards the 
vertical direction. We introduced the turbulence velocity which is equal to the sound velocity. 
However, we obtained only a small enhancement in the Ha equivalent width. Introduction of 
high velocity gradient may improve the situation. We calculated the non-LTE state for the 
model proposed by Stee & Araiijo (1994), keeping the disk mass of our standard disk model. 
However, the resultant Ha equivalent width does decrease drastically, since the high density 
region becomes compact. Sigut & Jones (2007) adopted the N{zu,0) oc law and obtained 

the Ha luminosity comparable to the observed values. We also made the preliminary calcu- 
lation for the N{w, 0) oc t<7~^'^ law with the same mass loss rate as in model 7 and obtained 
the Ha emission equivalent width of 80 A for i = 90°, reflecting the large projected area of 
the optically thick region at Ha. This result may suggest that the density in actual Be stars 
decreases more slowly than our model with A^(ci7,0) oc as suggested by Sigut & Jones 

(2007). This is possible if the disk temperature decreases outwards, caused by deficient Balmer 
continuous radiation in the optically thick region. 

Our adoption of the single-fiight escape probability for the lines may not be valid in the 
solution of the line radiative transfer. It is evident that this approximation underestimates the 
net radiative bracket, because the local escape of photons is replaced by the global escape from 
the disk. It is also important to get the exact solution numerically for the line transfer problem. 
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